Quantum Monte Carlo simulations of fidelity at magnetic quantum phase transitions 
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When a system undergoes a quantum phase transition, the ground-state wave-function shows a 
change of nature, which can be monitored using the fidelity concept. We introduce two Quantum 
Monte Carlo schemes that allow the computation of fidelity and its susceptibility for large interacting 
many-body systems. These methods are illustrated on a two-dimensional Heisenberg model, where 
fidelity estimators show marked behaviours at two successive quantum phase transitions. We also 
develop a scaling theory which relates the divergence of the fidelity susceptibility to the critical 
exponent of the correlation length. A good agreement is found with the numerical results. 
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What happens to the ground-state (GS) wave-function 
when a physical system goes across a quantum phase 
transition (QPT)? Rooted in quantum information the- 
ory, the fidelity approach [TJ [2] provides an interesting 
global answer in terms of the overlap between GS of the 
system at two different values of the driving parameter. 
The basic idea, which precursor may be found in Ander- 
son's orthogonality catastrophe [3] , is that close quantum 
states become more orthogonal close to a phase transi- 
tion. The resulting fidelity drop can then provide a useful 
probe to detect the QPT. This is particularly interesting 
as the fidelity is a global, model-independent quantity 
that incorporates all the information contained in the 
GS wave-functions. This is opposite to other usual ap- 
proaches to phase transitions, which often need an input 
such as the knowledge of a specific order parameter. 

Consider the Hamiltonian 

H = H + \H\ (1) 

with H\ acting as a perturbation to Hq. The fidelity 
-FXAj^Aa) is defined as the modulus of the overlap be- 
tween GS of H at two different values of A: 

F(A 1 ,A 2 ) = |(^|^ A2 )|. 

Suppose that the system undergoes a QPT for a value A c 
of the driving parameter. We expect the fidelity to have 
a singular behaviour when either Ai or A2 are close to A c , 
especially if the difference SX — A2 — Ai is small [1 . When 
<5 A — > 0, the fidelity is dominated by its leading term, the 
fidelity susceptibility xf 0], with F ~ 1 — ^~Xf- 

When H describes a many-body system, computations 
of F or xf are complicated. Besides a few analytical 
results on specific models [2], the main effort has been 
put in their numerical evaluation. Exact diagonalization 
(ED) and tensor-network (TN) methods 5 -including 
density matrix renormalization group (DMRG) [6 - have 
been the most widely used techniques in that respect, 
even though they suffer from several caveats. The ED 
method needs the full computation of the GS wave- 
function and is therefore limited to small systems. TN 



methods provide a variational ansatz for the GS wave- 
functions, allowing a straightforward computation of 
overlaps. This ansatz turns out to be excellent for one- 
dimensional systems, where DMRG [6] in particular has 
proved its full strength. Recently, several TN based 
works studied two-dimensional (2d) systems and their 
fidelity properties [7]. However, these methods remain 
variational and may fail in correctly capturing the GS 
properties of complex many-body Hamiltonians in d > 1 , 
especially close to a QPT. 

In this Letter, we present two different Quantum 
Monte Carlo (QMC) schemes which allow an exact (al- 
beit stochastic) computation of the fidelity F and its sus- 
ceptibility xf- We apply these methods to the antifcr- 
romagnetic (AF) Heisenberg spin model on a 2d lattice. 
Varying one exchange coupling in the spin model causes 
two successive phase transitions, both of which are found 
to be captured by the fidelity and its susceptibility. In 
passing, we derive a scaling theory for the divergence 
of xf at a second-order QPT. The two schemes benefit 
from the power of QMC methods, which allow to treat 
very large systems in any dimension. The first scheme, 
which calculates F, is applicable to all AF systems ad- 
mitting a singlet GS. The second scheme for xf is even 
more general and can be applied to several many-body 
problems. This opens a novel path for the fidelity ap- 
proach to QPT. These techniques are only efficient when 
the underlying QMC method is, i.e. when there is no 
sign problem. 

Model — The schemes are illustrated on the spin- 1/2 
Heisenberg Hamiltonian 

H = H + XH X = Y, s i ■ S J + A Si ' S i ( 2 ) 

{ij)o 

on the CaVO lattice [SJ, a l/5th depleted square lattice 
(see Fig. [I]). The first sum runs over nearest-neighbor 
spins on 0-bonds (solid lines in Fig. [T]) while the second 
is over A-bonds (dashed lines) . This lattice structure can 
be found in the AF compound CaV^g (hence the lattice 
name) , even though the interactions are more complex in 
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FIG. 1: (color online) (a) CaVO lattice with nearest neighbor 
bonds of type A (dashed lines) and (solid lines), (b)-(c) 
Typical VB states on the CaVO lattice, (d) Overlap graph of 
previous VB states forming two loops. 



the real compound [5j . Varying the coupling A allows the 
occurrence of two QPT separating an intermediate Neel- 
ordered AF phase from respectively a low-A plaquette 
and a high-A dimer phase [5] [TO] . H conserves the total 
spin of the system and in particular, AF interactions A > 
lead to a singlet GS. 

Fidelity measurements — Measuring the fidelity seems 
at first glance easy within a standard Projection QMC 
scheme [TJJ. Decomposing the ground-state |V>o) = 
Y^i a i\<Pi) m the simulation basis {| <£>)}, one generates 
representatives \ip^) {i.e. in proportion of |<ij|) of the 
GS via the projection scheme for two different values of 
A. The problem comes from the fact that the QMC es- 
timate of the fidelity (tp^ 1 \f^ 2 ) will vanish most of the 
time in the commonly used orthogonal basis, leading to 
a serious statistical problem. However, when the GS is a 
singlet, it can be decomposed in the Valence Bond (VB) 
basis, which has a crucial non-orthogonality property. In- 
deed, any two VB states always have a non-zero fidelity 
F = \{<pi\ip2)\ = 2 Ne ~ N / 2 where N e is the number of 
loops obtained by superimposing the two VB states (see 
Fig. [I]), and N the total number of spins. This prop- 
erty solves the statistical problem and allows an efficient 
computation of the fidelity. 

More specifically, we work with a VB projector loop 
algorithm recently proposed by Sandvik and Evertz |12j . 
To avoid the sign problem, we simulate non-frustrated 
AF on bipartite lattices, which leads to real positive val- 
ues of coefficients a* > and of VB overlaps (<pi\tp2) > 0. 
In the VB loop algorithm [T3], two VB representatives 
\<Pl) and \<px) of the ground-state are generated by prop- 
agating two initial VB states. Simulating at the same 
time two different physical systems with couplings Ai and 
A2 allows a QMC estimator of the square of the fidelity: 



^ 2 (Ai,A 2 ) 
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F(\i,\2) can be computed for any value of Ai and 
A2 for all models that can be simulated with VB QMC 
methods. 
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FIG. 2: (color online) Fidelity per site / as function of Ai and 
A2 for different system sizes. Expansion power n of H in the 
VB QMC method is n/N = 20 for N < 12 2 , 10 for N = 16 2 
and 4 for N > 24 2 . A range is [0, 2.5] for N < 16 2 and [0.5, 2] 
for N > 24 2 . Resolution AA for the plots is 0.02 for TV = 16 2 , 
0.01 otherwise. 



In the following, we illustrate this method for the 
Heisenberg model on the CaVO lattice (Eq. [2]) . The unit 
cell contains 4 spins, and we simulated square samples 
with L x L unit cells (total number of spins N = 4L 2 ) 
up to L — 16, using periodic boundary conditions. For 
such large systems, the fidelity essentially vanishes for 
all Ai 7^ A 2 . As suggested in Ref. [13] . we compute the 
fidelity per site /(Ai,A 2 ) = F(Ai, A 2 ) 1 ^ Ar , which is well- 
behaved as N — > 00. 

Our data for the fidelity per site are presented on 
Fig. [2] Around the diagonal where /(A, A) = 1, we no- 
tice the appearance of two pinch points, roughly around 
Aj G [0.8, 1.1] and A^ e [1.5, 1.8]. It has been argued that 
these features are characteristic of continuous QPT [T3] 
and our results are in agreement with the two well-known 
second-order QPT in this model. Far enough away from 
these two critical regions we notice that there is no sig- 
nificant change of / when N — * 00. Within the critical 
regions / drops faster with system size. Given our sta- 
tistical errors (up to 2% for the chosen range of A), we 
cannot however provide more precise ranges for the crit- 
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ical points. To locate more accurately the QPT, we now 
turn to the leading correction of fidelity around the di- 
agonal Ai = A 2 . 

Fidelity susceptibility - For SX — > 0, we con- 
sider the fidelity susceptibility xf which can be ex- 
pressed jl] as the imaginary-time integral xf — 
J °° t[(H x (0)H x (t)) - (#a(0)) 2 ] dr. This definition of- 
fers a natural extension to finite temperature T — 1//3: 

X f(P) = / t [(H x (0)H x (r)) - (H x (0)) 2 ] dr. (3) 
Jo 

and xf = hm/3^00 Xf(P)- This definition of xf{P) dif- 
fers from the Bures metric ds 2 usually defined for mixed 
states [T3], even though both have the same T = limit. 
However, one can prove [TS] that ds 2 /2 < xf(P) < ds 2 , 
showing that both quantities scale in the same way. 

An advantage of Eq. ^ is that xf(P) can be com- 
puted within a QMC stochastic series expansion (SSE) 
formalism [TBI EZ] ■ Note the importance of taking J3/2 as 
the upper limit of the integral as the /^-periodicity of the 
path integral would lead otherwise to incorrect results. 
In the SSE formalism the partition function is expanded 
in powers of /3, Z = X^^Lo TriJ". Starting from 

the expression of time-displaced correlations functions in 
SSE (Eq. 3.15 of Ref.QU, xf(P) is estimated as 



Xf(P) 



1 

A2 



J2 A(m,n)(N x (m)) - (N x 



(4) 



where N x (m) is the number of times two elements 
of H x appear separated by m positions in the SSE 
sequence |16j and the total number of appear- 
ance of elements of H x . The amplitude A(m,n) — 

(n-m-ajim! ^ d T T m+1 {I - r )n-n»-2 can be ca l cu l ate d 

for all (to, n) prior to simulations by numerical integra- 
tion or analytically for large n [15]. We emphasize that 
this formalism allows to compute Xf{P) f° r an Y model 
which can be simulated with SSE. 

The computation of xf(P) can turn costly for large 
systems at low T. We reached L < 16 and used (3 = 
10L for the CaVO lattice, and limited simulations to the 
relevant A range for the largest L and j3. 

Fig. [3ja)-(f) display the susceptibility fidelity per site 
Xf/N, showing the apparition of two peaks as a func- 
tion of A. While the N = 4 2 and 8 2 samples show rather 
broad feature (especially for the second peak), the peaks 
are clearly emergent as system size is increased and tem- 
perature lowered. From the position of the two peaks 
for the lowest T and largest size, one obtains estimates 
Aj = 0.94(1) and A 2 = 1.65(1) for the two quantum crit- 
ical points, in full agreement with QMC computations of 
order parameter and spin gap [10]. Note that the posi- 
tions of the maxima of Xf at finite T (see Fig. [3j^e)-(f ) ) 
also allow to determine faithfully A \ and A 2 , even though 
a small shift is observed if /3 is too low. We therefore find 
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FIG. 3: (color online) (a)-(f) Fidelity susceptibility per site 
Xf/N versus A for different system sizes iV and inverse tem- 
perature (3. (g) Xf/N versus A for the largest f3 for different 
TV. (h) Scaling of the two peaks of Xf/N versus linear size L 
(log-log scale). Lines denote power-law fits. 



that the fidelity susceptibility behaves as a good global 
indicator of QPTs in a 2d quantum system. 

Away from criticality, xf is extensive in all phases - 
see Fig. |3^g) . Scaling of the peaks in Fig. |3^h) reveals a 
power-law divergence at criticality Xi?(A c )/iV ~ U° , with 
lu = 0.73(3) for the first QPT and lu = 0.79(6) for the 
second (error bars originate from the statistical error bar 
in the QMC data). This agrees with the behaviour of / 
which develops pinch points at criticality but essentially 
does not change in non-critical regions as N — > 00. The 
observed symmetry of the peaks around their divergence 
explains the hourglass shape of the pinch points in /. We 
also understand why the second pinch point is harder to 
see in Fig. [2J indeed xf is much smaller close to A 2 than 
to \\ (almost a factor of 2 in all cases). 

Scaling theory — We now account for the divergence 
at the QPT by formulating a finite-size scaling (FSS) 
theory. FSS theories for xf have been proposed earlier 
but the relation to standard critical exponents at second- 
order QPT has been missed. In the most elaborate work, 
Campos Venuti and Zanardi [18J discuss the divergence of 
Xf as a function of the scaling dimension [H x ] of the part 
of the Hamiltonian that drives the transition, i.e. how 
this operator scales at A c : H x ~ L^^ Hx ^. Here we explic- 
itly calculate this scaling dimension. From the definition 
of the correlation length critical exponent £ ~ (A — A c ) _ly , 
we have [A] = \/v. Noting that [H] = z where z is the 
dynamical critical exponent, we deduce [H x ] = z — 1/u 
from Eq. Q. Finally, we conclude from Eq. ^ that 
[xf] = 2[t] + 2[H x ] = -2jv (see also Ref. [IS]) and there- 
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fore [xf/N] = —2/v + d, leading to the prediction: 

X F (\c)/N~L^- d . (5) 

This relation should hold for all second-order QPT and 
explains the superextensive behaviour generally observed 
for xf in terms of the usual critical exponents. From the 
value v ~ 0.7112 of the universality class of the 3d 0(3) 
model [19 to which belong both QPT studied here, we 
expect from this analysis lo ~ 0.812, in agreement with 
our numerical estimates. 

Discussions and conclusion — In conclusion, we pre- 
sented two QMC schemes that are able to calculate with 
high accuracy the fidelity and its susceptibility for quan- 
tum interacting systems, in any dimension. This allows 
to pin down the behaviour of fidelity at QPT, using one 
of the most sophisticated numerical techniques for the 
many-body problem. Taking the example of the Hciscn- 
berg model on the CaVO lattice, we find that both F 
and Xf are able to locate the two quantum critical points 
present in this system. The fidelity susceptibility acts as 
a more precise indicator as criticality manifests itself as a 
marked peak in xf- However, there is in principle more 
information contained in F. This could be useful to de- 
tect transitions that xf does not capture [2U] , as well as 
in the context of quantum quenches [211 EZj ■ 

The divergence of xf at criticality is accounted for 
by the scaling theory that we have presented, where the 
connection to the correlation length exponent of the uni- 
versality class of the QPT is made. We also showed that 
the generalization of xf to finite temperature (Eq. |3| al- 
lows to detect criticality for moderate values of T. This 
is of practical interest as simulations can be performed 
at a smaller computational cost. 

The method proposed for measuring xf works for any 
model which can be simulated within the generic SSE 
scheme [151 EZ] , opening the door to the study of fidelity 
in many different physical systems. We expect that our 
scheme can be extended to measure the Loschmidt echo, 
another witness of quantum criticality |23j , which can be 
measured experimentally in this context |24j . 

Note added — The scaling relation Eq. (|5| has been 
independently derived in recent preprints |22j . 
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